#params_det$s.mod.PFS <- "Generalised gamma" #input list is adapted
l.det.res <- CEModel_det(params = params_det) # run the function
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
#params_det$s.mod.PFS <- "Generalised gamma" #input list is adapted
l.det.res <- CEModel_det(params = params_det) # run the function
#params_det$s.mod.PFS <- "Generalised gamma" #input list is adapted
l.det.res <- CEModel_det(params = params_det) # run the function
params_det$s.mod.PFS <- "Gamma"
#params_det$s.mod.PFS <- "Gamma" #input list is adapted
l.det.res <- CEModel_det(params = params_det) # run the function
rm(list = ls())  # delete everything that is in R's memory
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
setwd(paste(path,"report", sep = "")) # set working directory to File Source
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
params_det$s.mod.PFS <- "Gamma"
save(ls(), file = paste(path,"data/01_model_inputs.R", sep = ""))
save(list = ls(), file = paste(path,"data/01_model_inputs.R", sep = ""))
rm(list = ls())  # delete everything that is in R's memory
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
rm(list = ls())  # delete everything that is in R's memory
options(scipen=1000) # removal of scientific notation
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
setwd(paste(path,"report", sep = "")) # set working directory to File Source
#load package
library(rms)
library(survival)
library(muhaz)
library(survminer)
library(ggplot2)
library(data.table)
library(printr)
library(dplyr)
library(flexsurv)
library(tidyverse)
library(scales)
library(matrixStats)
library(reshape2)
# load functions
source(paste(path,"R/01_fun_analyses.R", sep = ""))
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
params_det$s.mod.PFS
params_det$s.mod.PFS <- "Generalised gamma"
params_det$s.mod.PFS
#params_det$s.mod.PFS <- "Generalised gamma" #input list is adapted
l.det.res <- CEModel_det(params = params_det) # run the function
params_det$s.mod.PFS <- "Generalisedgamma"
#params_det$s.mod.PFS <- "Generalisedgamma" #input list is adapted
l.det.res <- CEModel_det(params = params_det) # run the function
Time <- c(0:abs(n.tw))*n.days.week/(365.25/12)
l.det.res$Trace.trt <- cbind(l.det.res$Trace.trt, Time)
l.det.res$Trace.comp <- cbind(l.det.res$Trace.comp, Time)
plot(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"PFS"],
type = "l",
lwd=2, col="blue", xlab = "Time (in months)",
ylab = "Propotion of patients in each state",
ylim = c(0,1),
cex.lab=0.75, main = "Eribulin")
lines(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"PD"],
type = "l",
lwd=2, col="red")
lines(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"Dead"],
type = "l",
lwd=2, col="black")
legend("topright", legend = c("PFS", "PD", "Dead"), col=c("blue", "red", "black"),
lty=c(1,1,1), cex=0.8)
plot(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"PFS"],
type = "l",
lwd=2, col="blue", xlab = "Time (in months)",
ylab = "Propotion of patients in each state",
ylim = c(0,1),
cex.lab=0.75,  main = "Non-eribulin")
lines(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"PD"],
type = "l",
lwd=2, col="red")
lines(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"Dead"],
type = "l",
lwd=2, col="black")
legend("topright", legend = c("PFS", "PD", "Dead"), col=c("blue", "red", "black"),
lty=c(1,1,1), cex=0.8)
params_det$s.mod.PFS <- "Gamma"
#params_det$s.mod.PFS <- "Gamma" #input list is adapted
l.det.res <- CEModel_det(params = params_det) # run the function
Time <- c(0:abs(n.tw))*n.days.week/(365.25/12)
l.det.res$Trace.trt <- cbind(l.det.res$Trace.trt, Time)
l.det.res$Trace.comp <- cbind(l.det.res$Trace.comp, Time)
plot(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"PFS"],
type = "l",
lwd=2, col="blue", xlab = "Time (in months)",
ylab = "Propotion of patients in each state",
ylim = c(0,1),
cex.lab=0.75, main = "Eribulin")
lines(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"PD"],
type = "l",
lwd=2, col="red")
lines(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"Dead"],
type = "l",
lwd=2, col="black")
legend("topright", legend = c("PFS", "PD", "Dead"), col=c("blue", "red", "black"),
lty=c(1,1,1), cex=0.8)
plot(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"PFS"],
type = "l",
lwd=2, col="blue", xlab = "Time (in months)",
ylab = "Propotion of patients in each state",
ylim = c(0,1),
cex.lab=0.75,  main = "Non-eribulin")
lines(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"PD"],
type = "l",
lwd=2, col="red")
lines(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"Dead"],
type = "l",
lwd=2, col="black")
legend("topright", legend = c("PFS", "PD", "Dead"), col=c("blue", "red", "black"),
lty=c(1,1,1), cex=0.8)
params_det$s.mod.PFS <- "Weibull"
#params_det$s.mod.PFS <- "Weibull" #input list is adapted
l.det.res <- CEModel_det(params = params_det) # run the function
Time <- c(0:abs(n.tw))*n.days.week/(365.25/12)
l.det.res$Trace.trt <- cbind(l.det.res$Trace.trt, Time)
l.det.res$Trace.comp <- cbind(l.det.res$Trace.comp, Time)
plot(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"PFS"],
type = "l",
lwd=2, col="blue", xlab = "Time (in months)",
ylab = "Propotion of patients in each state",
ylim = c(0,1),
cex.lab=0.75, main = "Eribulin")
lines(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"PD"],
type = "l",
lwd=2, col="red")
lines(x = l.det.res$Trace.trt[,"Time"],
y = l.det.res$Trace.trt[,"Dead"],
type = "l",
lwd=2, col="black")
legend("topright", legend = c("PFS", "PD", "Dead"), col=c("blue", "red", "black"),
lty=c(1,1,1), cex=0.8)
plot(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"PFS"],
type = "l",
lwd=2, col="blue", xlab = "Time (in months)",
ylab = "Propotion of patients in each state",
ylim = c(0,1),
cex.lab=0.75,  main = "Non-eribulin")
lines(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"PD"],
type = "l",
lwd=2, col="red")
lines(x = l.det.res$Trace.comp[,"Time"],
y = l.det.res$Trace.comp[,"Dead"],
type = "l",
lwd=2, col="black")
legend("topright", legend = c("PFS", "PD", "Dead"), col=c("blue", "red", "black"),
lty=c(1,1,1), cex=0.8)
rm(list = ls())  # delete everything that is in R's memory
options(scipen=1000) # removal of scientific notation
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
setwd(paste(path,"report", sep = "")) # set working directory to File Source
#load package
library(rms)
library(survival)
library(muhaz)
library(survminer)
library(ggplot2)
library(data.table)
library(printr)
library(dplyr)
library(flexsurv)
library(tidyverse)
library(scales)
library(matrixStats)
library(reshape2)
library(kableExtra)
# load functions
source(paste(path,"R/01_fun_analyses.R", sep = ""))
source(paste(path,"R/03_fun_dowsa.R", sep = ""))
source(paste(path,"R/03_fun_scenarios.R", sep = ""))
source(paste(path,"R/03_fun_psa_inputs.R", sep = ""))
source(paste(path,"R/03_fun_threshold.R", sep = ""))
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
df.scen <- as.data.frame(generate_scenario_results(params_det = params_det))
names(df.scen) <- c("Total QALY eribulin", "Total QALY non-eribulin", "Total costs eribulin", "Total costs non-eribulin", "Incremental Costs","Incremental QALY","ICER")
rownames(df.scen) <- c("Deterministic base-case", "Hospital perspective",
"Use of Weibull distribution to model OS",
"Use of log-logistic distribution to model PFS",
"Use of lognormal distribution to model TTD",
"Vial sharing for eribulin",
"Alternative utility values from ERG",
"No last year of life treatment",
"Use of hazard ratio from EMBRACE for OS and PFS")
df.scen[,c("Total costs eribulin", "Total costs non-eribulin", "Incremental Costs","ICER")] <- round(df.scen[,c("Total costs eribulin", "Total costs non-eribulin", "Incremental Costs","ICER")],0)
df.scen[,c("Total QALY eribulin", "Total QALY non-eribulin")] <- round(df.scen[,c("Total QALY eribulin", "Total QALY non-eribulin")],3)
df.scen <- df.scen[,c("Total QALY eribulin", "Total costs eribulin", "Total QALY non-eribulin", "Total costs non-eribulin", "Incremental QALY","Incremental Costs","ICER")]
knitr::kable(
df.scen, booktabs = TRUE,
format = "latex",
caption = 'Scenario analyses results') %>%
kable_styling(latex_options = "scale_down")
df.scen
knitr::kable(
df.scen, booktabs = TRUE,
format = "latex",
caption = 'Scenario analyses results') %>%
kable_styling(latex_options = "scale_down")
library(knitr)
l.det.res <- CEModel_det(params = params_det) # run the function
df.und.det.res <- l.det.res[[2]] # export results from the list to a separate dataframe
df.und.det.res[grep(pattern = "TC", x = rownames(df.und.det.res)),] <- paste('\u20ac',format(round(df.und.det.res[grep(pattern = "TC", x = rownames(df.und.det.res)),],0), big.mark = ",", scientific = FALSE)) # round costs to 0 decimal
df.und.det.res[grep(pattern = "TC", x = rownames(df.und.det.res), invert = TRUE),] <- round(as.numeric(as.character(df.und.det.res[grep(pattern = "TC", x = rownames(df.und.det.res), invert = TRUE),])),3) # convert LY and QALY to numbers & round LY and QALY to 3 decimals
rownames(df.und.det.res) <- gsub("TLY.", "Total life years - ", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub("TQ.", "Total QALY - ", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub("TC.", "Total costs - ", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub("pfs", " progression-free", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub("pd", " progressed disease", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".disutil", " disutility AE", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".systhx", " systemic treatment costs", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".adm", " administration", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".mon", " monitoring", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".surg", " surgery", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".hosp", " hospitalisation", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".radio", " radiotherapy", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".unrel.comp", " unrealted health care", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".travel", " travel", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub(".lastyear", " last year of life", rownames(df.und.det.res))
rownames(df.und.det.res) <- gsub("AE", "adverse event", rownames(df.und.det.res))
knitr::kable(
df.und.det.res, booktabs = TRUE,
caption = 'Undiscounted deterministic results'
)
df.scen <- as.data.frame(generate_scenario_results(params_det = params_det))
names(df.scen) <- c("Total QALY eribulin", "Total QALY non-eribulin", "Total costs eribulin", "Total costs non-eribulin", "Incremental Costs","Incremental QALY","ICER")
rownames(df.scen) <- c("Deterministic base-case", "Hospital perspective",
"Use of Weibull distribution to model OS",
"Use of log-logistic distribution to model PFS",
"Use of lognormal distribution to model TTD",
"Vial sharing for eribulin",
"Alternative utility values from ERG",
"No last year of life treatment",
"Use of hazard ratio from EMBRACE for OS and PFS")
df.scen[,c("Total costs eribulin", "Total costs non-eribulin", "Incremental Costs","ICER")] <- round(df.scen[,c("Total costs eribulin", "Total costs non-eribulin", "Incremental Costs","ICER")],0)
df.scen[,c("Total QALY eribulin", "Total QALY non-eribulin")] <- round(df.scen[,c("Total QALY eribulin", "Total QALY non-eribulin")],3)
df.scen <- df.scen[,c("Total QALY eribulin", "Total costs eribulin", "Total QALY non-eribulin", "Total costs non-eribulin", "Incremental QALY","Incremental Costs","ICER")]
knitr::kable(
df.scen, booktabs = TRUE,
format = "latex",
caption = 'Scenario analyses results') %>%
kable_styling(latex_options = "scale_down")
knitr::kable(
df.scen, booktabs = TRUE,
format = "latex",
caption = 'Scenario analyses results')
knitr::kable(
df.scen, booktabs = TRUE,
caption = 'Scenario analyses results')
tiff(file = paste(path,"figs/02_CEplane_paper.tif", sep = ""), res = 300, units= "cm", width = 10, height = 10)
ggplot(dfprob, aes(x = Inc.QALY, y = Inc.Costs)) +
geom_point(colour = "gray", shape = 1) +
geom_point(aes(x = mean(Inc.QALY), y = mean(Inc.Costs)), colour = "black", shape = 18) +
xlab ("Incremental QALY") +
ylab("Incremental costs") +
geom_hline(yintercept = 0, linetype="dashed",
color = "black") +
geom_vline(xintercept = 0, linetype="dashed",
color = "black") +
geom_abline(intercept = 0, slope = 80000, linetype= "dotted",
color = "black") + # 80,000 per QALY threshold line
stat_ellipse(type = "t",#assumes multivariate t distribution, using "norm would assume a multivariate normal distribution, 95%CI
colour = "black",
level = 0.95, size = 1) +
xlim(c(-limQALY,limQALY)) +
scale_y_continuous(labels = dollar_format(prefix = "\u20ac ", suffix = ""),
limits = c(-limCosts,limCosts)) +
theme_bw()
dev.off()
#determine path, where the model is located
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
tiff(file = paste(path,"figs/02_CEplane_paper.tif", sep = ""), res = 300, units= "cm", width = 10, height = 10)
ggplot(dfprob, aes(x = Inc.QALY, y = Inc.Costs)) +
geom_point(colour = "gray", shape = 1) +
geom_point(aes(x = mean(Inc.QALY), y = mean(Inc.Costs)), colour = "black", shape = 18) +
xlab ("Incremental QALY") +
ylab("Incremental costs") +
geom_hline(yintercept = 0, linetype="dashed",
color = "black") +
geom_vline(xintercept = 0, linetype="dashed",
color = "black") +
geom_abline(intercept = 0, slope = 80000, linetype= "dotted",
color = "black") + # 80,000 per QALY threshold line
stat_ellipse(type = "t",#assumes multivariate t distribution, using "norm would assume a multivariate normal distribution, 95%CI
colour = "black",
level = 0.95, size = 1) +
xlim(c(-limQALY,limQALY)) +
scale_y_continuous(labels = dollar_format(prefix = "\u20ac ", suffix = ""),
limits = c(-limCosts,limCosts)) +
theme_bw()
dev.off()
dfprob <- read.csv(paste(path,"output/03_probabilistic_results.csv", sep = ""), sep = ",")
limQALY  <- max(abs(dfprob$Inc.QALY)) #define limit Incremental QALY axis
limCosts <- max(abs(dfprob$Inc.Costs))#define limit Incremental QALY axis
tiff(file = paste(path,"figs/02_CEplane_paper.tif", sep = ""), res = 300, units= "cm", width = 10, height = 10)
ggplot(dfprob, aes(x = Inc.QALY, y = Inc.Costs)) +
geom_point(colour = "gray", shape = 1) +
geom_point(aes(x = mean(Inc.QALY), y = mean(Inc.Costs)), colour = "black", shape = 18) +
xlab ("Incremental QALY") +
ylab("Incremental costs") +
geom_hline(yintercept = 0, linetype="dashed",
color = "black") +
geom_vline(xintercept = 0, linetype="dashed",
color = "black") +
geom_abline(intercept = 0, slope = 80000, linetype= "dotted",
color = "black") + # 80,000 per QALY threshold line
stat_ellipse(type = "t",#assumes multivariate t distribution, using "norm would assume a multivariate normal distribution, 95%CI
colour = "black",
level = 0.95, size = 1) +
xlim(c(-limQALY,limQALY)) +
scale_y_continuous(labels = dollar_format(prefix = "\u20ac ", suffix = ""),
limits = c(-limCosts,limCosts)) +
theme_bw()
dev.off()
rm(list = ls())  # delete everything that is in R's memory
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
setwd(paste(path,"report", sep = "")) # set working directory to File Source
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
rm(list = ls())  # delete everything that is in R's memory
options(scipen=1000) # removal of scientific notation
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
setwd(paste(path,"report", sep = "")) # set working directory to File Source
#load package
library(rms)
library(survival)
library(muhaz)
library(survminer)
library(ggplot2)
library(data.table)
library(printr)
library(dplyr)
library(flexsurv)
library(tidyverse)
library(scales)
library(matrixStats)
library(reshape2)
library(kableExtra)
library(knitr)
# load functions
source(paste(path,"R/01_fun_analyses.R", sep = ""))
source(paste(path,"R/03_fun_dowsa.R", sep = ""))
source(paste(path,"R/03_fun_scenarios.R", sep = ""))
source(paste(path,"R/03_fun_psa_inputs.R", sep = ""))
source(paste(path,"R/03_fun_threshold.R", sep = ""))
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
View(df.pfs)
s.OS.trt  <- params_det$s.OS.trt
s.PFS.trt <- params_det$s.PFS.trt
s.TTD.trt <- params_det$s.TTD.trt
s.OS.comp  <- params_det$s.OS.comp
s.PFS.comp <- params_det$s.PFS.comp
s.TTD.comp <- params_det$s.TTD.comp
## extract s.mod.OS/PFS/TTD
s.mod.OS  <- params_det$s.mod.OS
s.mod.PFS <- params_det$s.mod.PFS
s.mod.TTD <- params_det$s.mod.TTD
time_pred <- seq(from = 0, to = 10, by = unit_time)
df.pfs$`L95%` <- as.numeric(as.character(df.pfs$`L95%`))
df.pfs$`U95%` <- as.numeric(as.character(df.pfs$`U95%`))
df.os$`L95%`  <- as.numeric(as.character(df.os$`L95%`))
df.os$`U95%`  <- as.numeric(as.character(df.os$`U95%`))
df.ttd$`L95%` <- as.numeric(as.character(df.ttd$`L95%`))
df.ttd$`U95%` <- as.numeric(as.character(df.ttd$`U95%`))
## extract s.OS/PFS and TTD estimation per treatment arm
s.OS.trt  <- params_det$s.OS.trt
s.PFS.trt <- params_det$s.PFS.trt
s.TTD.trt <- params_det$s.TTD.trt
s.OS.comp  <- params_det$s.OS.comp
s.PFS.comp <- params_det$s.PFS.comp
s.TTD.comp <- params_det$s.TTD.comp
## extract s.mod.OS/PFS/TTD
s.mod.OS  <- params_det$s.mod.OS
s.mod.PFS <- params_det$s.mod.PFS
s.mod.TTD <- params_det$s.mod.TTD
s.OS.trt[, "shape_Low"] <- 1 - pgamma(time_pred, shape = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape"), "est"]) * exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape(group1)"), "L95%"]), rate = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "rate"), "est"]) * exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "group1"), "est"]))
s.OS.trt[, "rate_Low"]  <- 1 - pgamma(time_pred, shape = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape"), "est"]) * exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape(group1)"), "est"]), rate = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "rate"), "est"]) * exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "group1"), "L95%"]))
s.OS.trt[, "shape_Up"]  <- 1 - pgamma(time_pred, shape = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape"), "est"]) * exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape(group1)"), "U95%"]), rate = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "rate"), "est"]) * exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "group1"), "est"]))
s.OS.trt[, "rate_Up"]   <- 1 - pgamma(time_pred, shape = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape"), "est"]) * exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape(group1)"), "est"]), rate = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "rate"), "est"]) * exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "group1"), "U95%"]))
s.OS.comp[, "shape_Low"] <- 1 - pgamma(time_pred, shape = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape"), "L95%"]), rate = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "rate"), "est"]))
s.OS.comp[, "rate_Low"]  <- 1 - pgamma(time_pred, shape = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape"), "est"]), rate = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "rate"), "L95%"]))
s.OS.comp[, "shape_Up"]  <- 1 - pgamma(time_pred, shape = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape"), "U95%"]), rate = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "rate"), "est"]))
s.OS.comp[, "rate_Up"]   <- 1 - pgamma(time_pred, shape = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "shape"), "est"]), rate = exp(df.os[which(df.os$Distnames == "Gamma" & df.os$Parnames == "rate"), "U95%"]))
s.PFS.trt[, "shape_Low"] <- 1 - pgamma(time_pred, shape = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape"), "est"]) * exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape(group1)"), "L95%"]), rate = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "rate"), "est"]) * exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "group1"), "est"]))
s.PFS.trt[, "rate_Low"]  <- 1 - pgamma(time_pred, shape = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape"), "est"]) * exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape(group1)"), "est"]), rate = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "rate"), "est"]) * exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "group1"), "L95%"]))
s.PFS.trt[, "shape_Up"]  <- 1 - pgamma(time_pred, shape = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape"), "est"]) * exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape(group1)"), "U95%"]), rate = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "rate"), "est"]) * exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "group1"), "est"]))
s.PFS.trt[, "rate_Up"]   <- 1 - pgamma(time_pred, shape = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape"), "est"]) * exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape(group1)"), "est"]), rate = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "rate"), "est"]) * exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "group1"), "U95%"]))
s.PFS.comp[, "shape_Low"] <- 1 - pgamma(time_pred, shape = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape"), "L95%"]), rate = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "rate"), "est"]))
s.PFS.comp[, "rate_Low"]  <- 1 - pgamma(time_pred, shape = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape"), "est"]), rate = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "rate"), "L95%"]))
s.PFS.comp[, "shape_Up"]  <- 1 - pgamma(time_pred, shape = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape"), "U95%"]), rate = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "rate"), "est"]))
s.PFS.comp[, "rate_Up"]   <- 1 - pgamma(time_pred, shape = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "shape"), "est"]), rate = exp(df.PFS[which(df.PFS$Distnames == "Gamma" & df.PFS$Parnames == "rate"), "U95%"]))
s.pfs.trt[, "shape_Low"] <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape(group1)"), "L95%"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "group1"), "est"]))
s.PFS.trt[, "shape_Low"] <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape(group1)"), "L95%"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "group1"), "est"]))
s.PFS.trt[, "rate_Low"]  <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape(group1)"), "est"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "group1"), "L95%"]))
s.PFS.trt[, "shape_Up"]  <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape(group1)"), "U95%"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "group1"), "est"]))
s.PFS.trt[, "rate_Up"]   <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape(group1)"), "est"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "est"]) * exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "group1"), "U95%"]))
s.PFS.comp[, "shape_Low"] <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "L95%"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "est"]))
s.PFS.comp[, "rate_Low"]  <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "est"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "L95%"]))
s.PFS.comp[, "shape_Up"]  <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "U95%"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "est"]))
s.PFS.comp[, "rate_Up"]   <- 1 - pgamma(time_pred, shape = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "shape"), "est"]), rate = exp(df.pfs[which(df.pfs$Distnames == "Gamma" & df.pfs$Parnames == "rate"), "U95%"]))
rm(list = ls())  # delete everything that is in R's memory
options(scipen=1000) # removal of scientific notation
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
setwd(paste(path,"report", sep = "")) # set working directory to File Source
#load package
library(rms)
library(survival)
library(muhaz)
library(survminer)
library(ggplot2)
library(data.table)
library(printr)
library(dplyr)
library(flexsurv)
library(tidyverse)
library(scales)
library(matrixStats)
library(reshape2)
library(kableExtra)
library(knitr)
# load functions
source(paste(path,"R/01_fun_analyses.R", sep = ""))
source(paste(path,"R/03_fun_dowsa.R", sep = ""))
source(paste(path,"R/03_fun_scenarios.R", sep = ""))
source(paste(path,"R/03_fun_psa_inputs.R", sep = ""))
source(paste(path,"R/03_fun_threshold.R", sep = ""))
# load model inputs
load(paste(path,"data/01_model_inputs.R", sep = ""))
#Draw tornado diagram
##SOURCE tornado diagram: https://stackoverflow.com/questions/55751978/tornado-both-sided-horizontal-bar-plot-in-r-with-chart-axes-crosses-at-a-given
df <- generate_dowsa(params_det = params_det)
source(paste(path,"R/03_fun_dowsa.R", sep = ""))
#Draw tornado diagram
##SOURCE tornado diagram: https://stackoverflow.com/questions/55751978/tornado-both-sided-horizontal-bar-plot-in-r-with-chart-axes-crosses-at-a-given
df <- generate_dowsa(params_det = params_det)
install.packages('dampack')
rm(list = ls())  # delete everything that is in R's memory
#determine path, where the model is located
path <- path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
source(paste(path,"R/01_fun_analyses.R"), sep = "")
load(paste(path,"data/01_model_inputs.R"), sep = "")
library(xlsx)
library(scales)
library(ggplot2)
dfprob <- read.csv(paste(path,"output/03_probabilistic_results.csv", sep = ""), sep = ",")
r.CEAC <- calcCEAC(Q.trt = dfprob$TQ_disc.trt, C.trt = dfprob$TC_disc.trt,
Q.comp = dfprob$TQ_disc.comp, C.comp = dfprob$TC_disc.comp, v.wtp = v.wtp)
dfCEAC <- as.data.frame(r.CEAC)
rm(list = ls())  # delete everything that is in R's memory
#determine path, where the model is located
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
source(paste(path,"R/01_fun_analyses.R"), sep = "")
#determine path, where the model is located
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
source(paste(path,"R/01_fun_analyses.R"), sep = "")
rm(list = ls())  # delete everything that is in R's memory
#determine path, where the model is located
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
setwd(paste(path,"report", sep = "")) # set working directory to File Source
source(paste(path,"R/01_fun_analyses.R"), sep = "")
rm(list = ls())  # delete everything that is in R's memory
#determine path, where the model is located
path <- "C:/Users/PouwelsXGLV/Documents/Maastricht/Eribulin CEA/Model final_Gamma/" #set location of model, to adapt!
setwd(paste(path,"report", sep = "")) # set working directory to File Source
source(paste(path,"R/01_fun_analyses.R", sep = ""))
load(paste(path,"data/01_model_inputs.R", sep = ""))
library(xlsx)
library(scales)
library(ggplot2)
dfprob <- read.csv(paste(path,"output/03_probabilistic_results.csv", sep = ""), sep = ",")
r.CEAC <- calcCEAC(Q.trt = dfprob$TQ_disc.trt, C.trt = dfprob$TC_disc.trt,
Q.comp = dfprob$TQ_disc.comp, C.comp = dfprob$TC_disc.comp, v.wtp = v.wtp)
dfCEAC <- as.data.frame(r.CEAC)
dfCEAC
install.packages("dampack")
library(devtools)
devtools::install_github("DARTH-git/dampack")
library(fansi)
remove.packages("fansi", lib="~/R/win-library/3.6")
devtools::install_github("DARTH-git/dampack")
library(dampack)
library("dampack")
install.packages(c("bookdown", "knitr", "tinytex"))
library(bookdown)
knitr::opts_chunk$set(message = FALSE)
library(bookdown)
library(tinytex)
